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Quasars have long been known to be variable sources at all wavelengths. Their optical vari¬ 
ability is stochastic, can be due to a variety of physical mechanisms, and is well-described 
statistically in terms of a damped random walk modePl The recent availability of large col¬ 
lections of astronomical time series of flux measurements (light curvesjP® offers new data 
sets for a systematic exploration of quasar variability. Here we report on the detection of 
a strong, smooth periodic signal in the optical variability of the quasar PG 1302-102 with a 
mean observed period of 1,884 ± 88 days. It was identified in a search for periodic variability 
in a data set of light curves for 247,000 known, spectroscopically confirmed quasars with a 
temporal baseline of ~ 9 years. While the interpretation of this phenomenon is still uncer¬ 
tain, the most plausible mechanisms involve a binary system of two supermassive black holes 
with a subparsec separation. Such systems are an expected consequence of galaxy mergers 
and can provide important constraints on models of galaxy formation and evolution. 

Subparsec supermassive black-hole (SMBH) binary systems are not resolvable except possi¬ 
bly with long baseline radio interferometry. An alternative approach to their detection is through a 
modulated variability - caused by, for example, perturbations in their accretion disks or precession 
of relativistic jets, if they are present (see Eig. 1). The best known candidate, OJ 28"^, has shown a 
pair of outburst peaks every 12.2 years for at least the past century: this object can be interpreted as 
a secondary SMBH perturbing the accretion disk of the primary SMBH at regular intervald^. Sys¬ 
tematic searches for equivalent systems to dateP^^have attempted to identify them from broad-line 
velocity offsets in their optical and near-infrared spectra but cannot detect the closest pairs (with 
<0.1 pc separation). 

We applied a novel joint wavelet and autocorrelation function (ACE) based technique that 
identifies objects exhibiting strongly periodic behavior in their light curves (M.J.G. et al, manuscript 
in preparation) to the largest set of quasar time series currently available. These are drawn from 
the Catalina Real-time Transient Survey (CRTSjnEl. PG 1302-102 (see Eig. 2) is the strongest 
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periodic candidate out of twenty objects meeting the selection criteria: strong constant wavelet 
peak, strong ACF detection of periodic behaviour, sufficient temporal coverage for 1.5 or greater 
cycles at the detected period, and a phased light curve well-described by a sinusoid. For statistical 
comparison, we have also generated a simulated light curve for each known quasar based on a 
damped random walk (DRW) model, a standard statistical description of the optical variability of 
quasarP, using the CRTS time sampling. We find that only one object from the simulated data of 
247,000 quasars satisfies the same selection criteria, showing that the number of quasars selected is 
statistically significant and that strongly periodic behavior is not expected as an artifact of a DRW 
process. 

PG 1302-102 has a median C-band magnitude of 15.0 and a redshift of 0.2784 which gives 
an absolute C-band magnitude of My = —25.81, assuming the 9-year WMAP cosmolog}^^. It 
is outside the footprint of the Sloan Digital Sky Survey but is associated with bright infrared and 
X-ray sources. It is also a very bright (720 mJy at 4.86 GHz), core-dominated flat spectrum radio 
source. Its optical/NIR spectrum (see Fig. 3) shows broad emission lines (H/3, Ha, Pa-/3, Pa-a) 
with an inferred mass of \og{M / Mq) = 8.3 — 9.4 and the object appears to be radiating at or close 
to its theoretical Eddington limit {\og{L/LEdd) = 0). The light curve for the quasar is well-fit by 
a sinusoid with an observed period of 1, 884 ± 88 days (corresponding to a restframe period of 
1,474 ± 69 days) and an amplitude of ~0.14 mag. CRTS data (covering ~1.8 cycles; that is. May 
2005 to present day) is augmented by archival monitoring datdi^IIl available back to May 1993 
giving a total of 4.1 cycles’ coverage. This data is consistent with the behaviour seen in the past 
nine years of CRTS data, particularly with stochastic photometric variation imposed on a periodic 
signal. Further simulations show that the detection is statistically significant, with an observed 
signal 40 times the scatter from the mean. 

As PG 1302-102 is bright and nearby, it has featured in a number of studies of quasars 
and their host galaxies. The radio and optical structure of the source is noted to be unusual. 
Hubble Space Telescope (HST) imaging shows that the quasar resides in a luminous elliptical 
host, as typical for radio-loud quasar^. There are also two companion galaxies that lie at projected 
distances of 3 and 6 kpc. Several features in radio images of PG 1302-102, such as the small 
radio core and sharp bends in the radio structure very close to the central source, correspond with 
features seen in the opticaP^. An interpretation is that the host galaxy is a fairly old merger, but 
that there might be more recent activity with the radio source just turning on and possible radio 
jets just emerging from the host galaxy. There may also be some indication of relativistic beaming 
connected with a jet. It should be noted that OJ 287 also exhibits a similar radio and optical 
morphologjJ^. 

PG 1302-102 was spectroscopically monitored over a six-month period in 1 99Cp2 and showed 
no detectable (greater than 5 - 10%) change in any component of its spectrum over that time. 
This lack of variation is not inconsistent with the ~60 month period that we have identified. A 
SMBH binary may also exhibit double-peak broad line profiles in its spectrum for a small window 
of separation between the paii^ (although disk emission from an accretion disk around a single 
source may also produce the same effecP^. At closer distances, the two black holes dynamically 
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affect the broad-line region clouds as a single complex entity producing single-peaked spectral 
lines with asymmetric line profiles. The Balmer and Paschen series spectral lines in PG 1302-102 
do not show a double peak profile but are consistently asymmetric (see Fig. 4). In particular, a 
small bump on the red wing of H/3 has been reported^, implying a velocity shift of order 200 
km s“^ between the narrow and broad components of H/3. One proposed explanation for this is a 
binary system. 

The physical interpretation of the periodicity is uncertain, although its sinusoidal nature sug¬ 
gests that it is kinematic in origin: we consider three possibilities here: (1) The optical flux could 
be the superposition of thermal emission from the accretion disk and a non-thermal contribution 
from a processing jet, and such a model can fit the observed data (see Methods). The expected 
precession period with a single SMBH is p ~ 10^'^ — 10®'^ year^, much longer than the ob¬ 
served period. Thus, a binary SMBH origin for the jet precession is more plausible. In this latter 
case, the jet could process either as a result of inner disk precession due to the tidal interaction of 
an inclined secondary SMBH or because of the precession of a circumbinary disk warped by the 
SMBH binary. (2) Another possibility is a temporary hot spot in the inner region of the accretion 
disk, but this leads to implausible SMBH mass estimates of \og{M/M q) = 11.4 — 12.2, depending 
on the degree of rotation of the SMBH (the largest reported SMBH massed are of the order of 
log(M/M0) 10). However, with a SMBH binary, periodic mass accretion rates can give rise to 
an overdense lump in the inner circumbinary accretion diskP^. The spectral energy distribution of 
a circumbinary disk also has a steeper power lavP^ and so accretion variations will have a more 
noticeable effect at shorter wavelengths. (3) Yet another possibility is a warped disk eclipsing part 
of the continuum as it processes, although SMBH binaries are proposed as a possible cause for 
themP^. We note as well that light curves for objects known to exhibit these phenomena do not 
resemble that of PG 1302-102 (see SOM). 

If PG 1302-102 were to be described as a binary SMBH pair with a total virial mass of 
log(M/Mo) ~ 8.5, the observed period gives an upper limit separation of ~0.01 pc between 
the pair. This would mean that the system has evolved well into the “final parsec” scale. The 
expectation is that most binary SMBH systems will spend the majority of their lifetime at such 
separations (0.01 - 1 pc), in an intermediate phase of evolution between scattering any stars in the 
nuclear region and gravitational radiation dominanceP^. 

Further observations could test the different interpretations mentioned above, particularly 
reverberation mapping to measure the behavior of emission line response to continuum variations, 
which is expected to be different for alternate explanation^. Continued monitoring by CRTS 
and other synoptic surveys will obviously track future cycles, and historical photometric data from 
photographic plate collections may provide more data for previous ones. With decadal baselines, 
the predicted change in period of the system may be detectable. Future spectroscopic observations 
can also test whether the line asymmetries vary on binary orbital timescales. Multiwavelength 
observations should provide more information about the innermost regions of the quasar and the 
nature of the jet. The relationship between PG 1302-102 and its two nearby companions may also 
furnish insight into the merger history of this source, particularly as these may contain similarly 
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sized SMBHs to the secondary. Finally, if PG 1302-102 is a SMBH binary, it is a strong candidate 
for any gravitational wave experiment sensitive to nanohertz frequency waves, such as those using 
pulsar timing arrays as well as any future space-home gravitational wave detection mission. 
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Figure 1 The parameter space of SMBH binary pairs. The expected orbital periods for 
SMBH close binary pairs at the specified separations as a function of total black-hole 
mass. The solid upper line for each separation indicates a z = 5 track and the solid lower 
line a z = 0.05 track whilst the two internal dotted lines show = 1.0 (lower) and = 2.0 
(upper) tracks respectively. The hatched region indicates the range over which CRTS has 
temporal coverage of 1.5 cycles or more of a periodic signal. The pink shaded region 
shows the region of detection for the best CRTS candidate given the range of virial black- 
hole masses reported in the literature. Also shown (solid black star) is the location of the 
best known SMBH binary candidate, OJ 287®. 
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Figure 2 The composite light curve for PG 1302-102 over a period of 7,338 days (~ 20 
years). The light curve combines data from two CRTS telescopes (CSS and MLS) with 
historical data from the LINEAR and ASAS surveys, and the literature (see Methods for 
details). The error bars represent one standard deviation errors on the photometry values. 
The dashed line indicates a sinusoid with period 1,884 days and amplitude 0.14 mag. 
The uncertainty in the measured period is 88 days. Note that this does not reflect the 
expected shape of the periodic waveform which will depend on the physical properties of 
the system. MJD, modified Julian day. 
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Figure 3 The composite spectrum for PG 1302-102. This combines an archival GALEX 
spectrum (ultraviolet) with optical/near-infrared spectra taken with the Keck and Palomar 
200” telescopes in April and June 2014. Fx, flux density. The prominent emission lines 
are indicated. The median flux errors are 5.6 x 10~^® erg s ^ cm“^ for the GALEX data 
(A < 0.3/im), 4.5 x 10'^® and 7.6 x 10^^’^ erg s~^ cm'^, respectively, for the blue (0.3/im < 
A < 0.5/im) and red (0.5/im < A < 0.9/im) optical spectra from Palomar, and 4.6 x 10~^® 
erg s~^ cm'^ for the Keck near-infrared (A > 0.9/im) spectrum. 
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Figure 4 The profiles of the Balmer and Paschen series lines of PG 1302-102. The data 
have been modelled with a multi-component line fitting technique (see Methods for de¬ 
tails). a, h Balmer H/3 (a) and Ha (b) have been fit using a narrow component (dashed 
blue line) and a broad Gaussian (solid orange line). The dashed green line shows the 
linear continuum component, and the total fitted profile is shown as a solid red line. H/3 
requires a single Gaussian offset from the narrow component but Ha requires two com¬ 
ponents - a central Gaussian plus a red wing. c,d The Paschen lines (Pa/3 (c) and Paa 
(d)) also show a consistent small asymmetry on the red side. 
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Methods 


1 CRTS data 

The Catalina Real-time Transient Survey (CRTS) makes use of the Catalina Sky Survey whieh 
began in 2004, operated by Lunar and Planetary Laboratory at University of Arizona, and uses 
three teleseopes (designated CSS, MLS and SSS) to eover the sky between deelination 6 = —75° 
and +65° (~ 80% of the sky) in order to diseover near-Earth objeets and potentially hazardous 
asteroids. The full Catalina surveys data set eontains time series for approximately 500 million 
sourees to a limiting magnitude of U ~ 20 with an average of ~ 250 observations over a baseline 
of 9 years per souree. CRTS operates an open data policy and the data are publicly available at: 
http ://catalinadata.org, 

For photometric calibration we combine observations taken with all three Catalina Sky Sur¬ 
vey telescopes since no difference was found between these systems. This is not surprising since 
each telescope specifically uses the same type of 4k x 4k CCD camera and observations with all 
three telescopes are calibrated using the same software pipeline. All Catalina observations are 
transformed to Johnson V based on 50-100 stars selected as G-type stars using 2MAS;^ colours. 
For bright stars, this photometry provides repeated photometry accurate to ~ O.Od^. However, as 
the photometry is unfiltered, there are significant variations with object colour. 

The 2007 Fandolt UBVRI standard star catalo^^ provides 109 stars centered near declination 
—50° in the magnitude range 10.4 < U < 15.5 and in the color index range —0.33 < {B — V) < 
1.66, while the 2009 Fandolt catalo^^ provides 202 standard stars along the celestial equator in the 
magnitude range 8.90 < U < 16.30, and the colour index range —0.35 < {B — V) < 2.30, along 
with 393 standard stars from previous standard star catalogs. In total there are 445 Catalina light 
curves matching Fandolt standards. On average each standard is measured 134 times. A handful 
of stars that appeared to exhibit significant variability were removed. 

Median magnitudes were calculated for each light curve and these were then used to deter¬ 
mine the following transformation equations between Johnson V and Catalina Voss'- 


V = Vcss + 0.31 X (B - V)2 + 0.04, 

(1) 

V = Vcss + 0.91 X (V - Rf + 0.04, 

(2) 

V = Vcss + 1.07 X (V - 1)2 + 0.04. 

(3) 


The dispersion in the fits to these transformations are 0.059, 0.056 and 0.063 magnitudes, 
respectively, for U < 16. 
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2 Candidate selection 


We applied the weighted wavelet transformP^l and the z-transform diserete correlation function 
(ZDCpj^to the CRTS light curves of spectroscopically confirmed quasars. Both of these algo¬ 
rithms can detect (quasi-)periodic behaviour in irregularly sampled data. We define the period of 
the quasar from the largest peak in the ZDCF between the second and third zero-crossings of a 
Gaussian process model fit to the ZDCI®. We have verified that this agrees with periods deter¬ 
mined by other methods. The period uncertainty is defined as: 

1.483 X MAD 

ap = - , — 

where MAD is the median of the absolute deviations from the median of the time intervals between 
successive peaks in the ZDCF and N is the total number of peaks considered. 

We only considered those objects whose wavelet peak significance places them in the top 
quartile of the data set (in terms of significance). Most kinematically-caused variations should 
manifest as a (near)-Keplerian signal and so we also use the r.m.s. scatter around an expected 
sinusoidal waveform (best-fit truncated Fourier series of up to 6 terms) with the ZDCF period to 
identify those objects most closely exhibiting the expected behavior. We excluded those quasars 
where the scatter is greater than the la lower limit on the median absolute deviation of the light 
curve, i.e., r.m.s./mad > 0.67, since we need to account for the intrinsic variability of the quasar 
as well. We also restricted our selection to candidates with temporal coverage of more than 1.5 
cycles, assuming the ZDCF period. 

Simulated data For each quasar in our data set, we have generated a simulated light curve assum¬ 
ing that it follows a damped random walk (DRW) model. Using the actual observation times, ti, 
we replace the observed magnitudes with those that would be expected under a DRW model. The 
magnitude X (t) at a given timestep At from a previous value X{t — At) is drawn from a Gaussian 
distribution with mean and variance?^ 


E{X{t)\X{t - At)) = - At) + 6r(l - 


Var(X(t)|A'(« - At)) = ^[1 - 

We add a Gaussian deviate normalized by the photometric error associated with the magnitude to 
be replaced at each time t to incorporate measurement uncertainties into the mock light curves. 
For each light curve, we set br to its median value and use the rest frame DRW fitting functiond^ 
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log / = ^ + B log + C(M + 23) + D log 

where (A, S, (7, D) = (—0.51, —0.479, 0.113, 0.18) for / = SFoo = ra"^ 12 and (A, 5, C, H) = 
(2.4, 0.17, 0.03, 0.21) for f = t. M is the absolute magnitude of the quasar and Xrf is the rest- 
frame wavelength of the filter. The mass of the blaek hole is either the measured virial masd^ or 
is drawn from a Gaussian distributioiP^ 


p(}ogMBH\M) 


!— exp 
y/^a 


(log Mbh -1^)^ 

2a2 


where /i = 2.0 — 0.27M and a = 0.58 + O.OllM. 

Statistical significance To assess the statistieal signifieanee of the deteetion of PG 1302-102, we 
generated 1,000 simulated light eurves as above. From these we determined the mean weighted 
wavelet power speetrum as a funetion of time and frequeney and its varianee. The dominant 
signal in the observed WWZ speetrum of PG 1302-102 is then seen to be ~ 40 times above the 
eorresponding mean DRW value in terms of the expeeted standard deviation. 

We have also performed a periodieity analysis of the light eurve of PG 1302-102 using the 
generalized Lomb-Seargle method whieh shows a statistieally signifieant peak at the same period 
identified by the wavelet and autoeorrelation analyses. The false alarm probability is << 10“^^ - 
the 10“^^ level is P{uj) = 0.335 and the observed peak is at P{uj) = 0.818. 

Theoretical predictions Simple disk models for eireumbinary gas and the binary-disk interaetion 
have been usecP^ to eonsider the number of SMBH binaries expeeted in a variety of surveys, 
assuming that sueh objeets are in the final gravitational wave-dominated phase of eoaleseenee (this 
equates to separations less than ~0.01 pe for a 10®Mq SMBH binary). This approaeh has been 
eombinecP^ with merger tree assembly models to similarly prediet the number of expeeted SMBH 
binaries at wider separations where speetral line shifts may be seen (this equates to separations 
greater than ~ 0.2 pe for a 10 ®Mq SMBH binary). The latter shows that in a sample of 10,000 
quasars at 2 ; < 0.7, there should be ~10 objeets and this number inereases by a faetor of ~ 5 — 10 
for 2 ; < 1. We note, however, that these theoretieal arguments are still subjeet to eonsiderable 
uneertainties; for example, if the final parsee problem eannot be resolved then there will not be any 
binaries in the ~0.01 pe regime. 

Assuming a limiting magnitude of 17 ~ 20, a deteetable range of orbital periods from 20-300 
weeks (spanning both GW- and gas-dominated regimes), a survey sky eoverage of In steradians, 
and a redshift range of 0.5 - 4.5, we would expeet 450 SMBH binaries following these approaehes. 
Our finding of 20 eandidates from a sample of 240,000 quasars is therefore eonservative. 89,000 
quasars in our sample also have virial blaek hole mass estimated^ (23% at z > 2) and if we as¬ 
sumed that eaeh of these was a SMBH binary with a separation of 0.01 pe then the CRTS temporal 
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baseline is suffieient to deteet 1.5 eyeles or more in 63% of them (ineluding 55% of the z > 2 pop¬ 
ulation). Our seareh is therefore sensitive to a large fraetion of the elose SMBH binary population. 

We note that our approaeh assumes that periodieity assoeiated with SMBH binaries manifests 
in a Keplerian form. If there is a larger set of non-Keplerian periodie SMBH binaries, either flaring, 
sueh as OJ 287, or not, then the 20 objeets we have identified may be a small sample of the total 
elose binary SMBH population. 

3 Archival data 

LINEAR dataP These were ealibrated with pre-release photometry from Pan-STARRS using the 
g, r, and i bandpasses. Comparison stars with instrumental magnitudes (eedjnag) between 14 
and 17 were seleeted within 0.1 degrees of PG 1302-102 in LINEAR images, g — i eolours were 
used to eompute an r-band eorreetion so that a ealibration star with g — i = 0 has an instrumental 
magnitude, eed mag = r. Zeropoints for eaeh frame were then derived based on these stars. The 
reported bandpass for the ealibrated magnitudes is therefore approximately r. Magnitude errors 
are eomputed by SEXTRACTOR, with typieal r.m.s. errors between frames of ~ 0.1. 

ASAS data The nominal limiting magnitude for ASAS’ is / ~ 13 and so PG 1302-102 is very 
elose to the deteetion threshold. The low signal-to-noise for sueh an objeet is the primary eause of 
the large degree of seatter seen in ASAS data for this souree. 

Historic data Such data for PG 1302 from previous quasar monitoring campaigns is available in 
the literaturell^ll^. To put all data on the same photometric scale, offsets were applied to account 
for differences in the photometric systems used. Region of temporal overlap between a pair of data 
sets were used to derive offsets so that both data had the same median value. Where no temporal 
overlap exists, the phased light curve was used to determine the median offset. 

Earlier individual photometric observations also exist of PG 1302-102 but the observational 
errors on these are typically ~ 0.1 mag and so it is difficult to determine whether they agree with 
the extrapolated behaviour. They also tend to be in different passbands which requires color terms 
to convert to the V -passband to which CRTS is calibrated. However, colour terms for quasars are 
known to vary (“bluer when brighter”) so a constant value cannot be assumed (the quoted {B — V) 
values for PG 1302-02 have a range of at least 0.2 mag), which introduces an additional error to 
the transformed magnitude. Such historical data are thus of limited utility. 

4 Spectroscopic data 

Optical An optical spectrum was obtained using the Double Spectrograph on the Hale 200-inch 
Telescope at Palomar Observatory on UT 2014 April 22. We obtained two 250 s exposures in 
cloudy conditions using the 1.0” wide slit, the 5500 A dichroic, the 600 grating on the 

blue arm (Ablaze = 4000 A), and the 316 grating on the red arm (Ablaze = 7500 A). 
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On UT 2014 May 26, we obtained additional speetroseopy of PG 1302-102 using the Low 
Resolution Imaging Speetrometer (LRi;^^ on the Keek I teleseope. We obtained two 300 s ex¬ 
posures in non-photometrie eonditions, using the 1.5” slit, the 5600 A diehroie, the 600 
grism on the blue arm (Ablaze = 4000 A) and the 400 ^mm~^ grating on the red arm (Ablaze = 
8500 A). 

These set-ups provided moderate resolution speetra aeross the entire optieal window, A A3100 A- 
1 /xm. The data from both teleseopes were redueed using standard proeedures and ealibrated using 
observations of standard stars obtained on the same (non-photometrie) nights. 

Near-infrared We obtained a NIR speetrum of PKS1302-102 with the TripleSpee instrument on 
the Palomar 200-ineh Hale teleseope on UT 2014 April 15. Conditions were elear and the seeing 
was ~1 aresee. The souree was observed at an airmass of 1.3890 and was observed for four 300- 
seeond exposures in an ABBA dither pattern for a total of 20 minutes of on-souree exposure. A 
speetrum of tellurie standard AOV star HD 112304 was obtained immediately following the souree 
speetrum at an airmass differenee of 0.15. We redueed the data using a modified version of the 
Spextool data reduetion paekageP^ whieh also performs a tellurie eorreetion using the standard star 
speetrurrl^. 

Fig. 3 shows the eombined optieal and near-infrared speetrum, in whieh we sealed the (non- 
photometrie) optieal speetrum by a faetor of 2.86 to meet the near-infrared speetrum in their mu¬ 
tually overlapping wavelength region (0.96 - 1.05 /xm). The ehange in flux is likely primarily due 
to a eombination of weather variations and slit losses. The Balmer lines, H/) and Ha, and Pasehen 
lines, Pa-/3 and Pa-a, are marked with vertieal dashed lines. 

5 System parameters 

To estimate the (total) blaek hole mass for this souree, we used the standard method of single 
epoeh virial BH estimation!^ and adapted relations derived for Pasehen lines in the near- infrared^. 

To determine the full-width at half maximum (FWHM) of the broad line eomponent of the four 
lines marked in Figure 4, we applied a multi-eomponent line fitting teehnique!531Zl in whieh we 
model the narrow-line eomponent of the line profiles by first fitting the [OIII] 4959,5007A lines 
and fixing the width of a narrow line eomponent in eaeh profile. For Ha we also inelude the 
[NII]6548,6583A doublet, fixed at a flux ratio of 2.96. The broad eomponent ean be modelled by 
up to three Gaussians. If the ratio of ehi-squares for sueeessive fits is greater than 0.8, an additional 
eomponent is added. The NIR speetrum has an error array whieh is used in the line modelling and 
parameter estimation (see below). For the optieal speetrum, an error array is estimated from the 
median differenee between adjaeent pixels (R. White, personal eommunieation). 

The measured values from the 2014 data are: 4,450 ± 150 km s“^ for H/), 2,520 ± 30 km s“^ 
for Pa-a, and 3, 200 ± 20 km s“^ for Pa-/). Errors on the FWHM were eomputed using a Monte 
Carlo approaeh: the best-fit model for a line was perturbed with a random draw from the error 
array at eaeh wavelength element and a new fit made. This was repeated 100 times for eaeh line 
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and the standard deviation of the broad line eomponent FWHM used as the error. These values give 
a (total) mass for the SMBH of log(M/M 0 ) = 8.8 ±0.6 (H/3), 8.5 ±0.1 (Pa-/?), and 8.4 ±0.1 (Pa- 
a). Previously reported values of the FWHM of H/? give a range of estimates for the SMBH mass 
in the literature using various teehniques of log(M/MQ) = 8.3 — 9.4 so our results are eonsistent 
with these. 

We note that the speetral fits of H/? and Pa-/? are not perfeet. However, the uneertainties in 
the FWHM and eontinuum luminosity that this introduees are small eompared to the broad range 
of SMBH mass estimated from the different lines. 

6 Alternate interpretations 

We present further diseussion here on the alternate interpretations eonsidered. We note that none 
of the objeets mentioned pass our eandidate seleetion eriteria. 

Jet related The optieal flux eould be the superposition of thermal eontribution from the aeeretion 
disk with non-thermal eontribution produeed by the underlying jet. The flux density of the jet (in 
the optieally thin regime) will be boosted in the observer’s frame relative to the eomoving frame: 


s,(v) = 

where a is the speetral index {Sj{u) oc z/“") and p = 2 for a eontinuous jet. The jet is preeessing 
with eonstant angular veloeity u, has an opening angle and an axis defined by the angles 00 
(between the jet axis and the line of sight) and tjq (the position angle in the plane of the sky): 


sin^ 0 = (sin cos ut ± cos sin 0o sin ± (sin cos 0o sin ut ± cos sin 0o cos poY 

Assuming a constant Lorentz factor 7 for the relativistic bulk motion of the jet, 7 = (1 — /?^)“^/^, 
and the Doppler factor is: 6 = 7“^(1 — /? cos 0)“^ Modeling the light curve in this way, we get 
best-fit parameters of: 7 = 5.4 ± 0.1, D = 0.5° ± 0.1, 0o = 5.0° ± 0.2, and rjo = 0.6° ± 1.4 
(assuming a = 1.66). 

A number of radio-loud quasars have been reportecP^® as showing periodic variability in 
their radio light curves. Whilst a SMBH binary could explain this, a more likely explanation is 
shock interaction with a helical jet or precession of a jet. However, the optical light curves of these 
objects (see Extended Data Fig. 1) do not show the distinctive behavior seen in that of PG 1302- 
102 suggesting that a different physical mechanism is more likely. We note as well that of the 20 
objects in our full sample showing optical periodicity, only 3 are associated with a radio source. 

Warped accretion disk Warped disks have been observed in a handful of AGhP®^ and the sug¬ 
gestion here is that as a warp precesses, it could obscure a small amount of continuum emission 
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which would then appear quite regular. Again there is no indieation of any periodic behavior in 
the CRTS light eurves available for known objeets with warped disks (see Extended Data Fig. 2) 
similar to that seen in PG 1302-102. 

PG 1302-102 shows a 14% variation in flux whieh would suggest that the size of the warp in 
the disk is quite large. This would also be an orientation-dependent phenomenon and as the souree 
is a blazar, its aeeretion disk should be oriented elose to faee-on to us and so any obseuring faetor 
should be limited in effeet. We also note that many stellar systems with warped aeeretion disks are 
resolvable binary systems. 
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Figure 5 Extended Figure 1 The optical light curves of quasars showing radio periodicity. 
Shown are the CRTS light curves for 11 quasars reportedSEHto show periodicity in their 
radio emission. Each light curve has been normalized to zero mean and individual curves 
are offset by a constant of 1.5 mag from each other. The data are split across two panels 
for ease of viewing. Error bars shown are standard la photometric errors. The CRTS light 
curve of PG 1302-102 (solid black stars) is also shown for comparison. 
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Figure 6 Extended Figure 2 The optical light curves of quasars with warped accretion 
disks. Shown are the CRTS light curves for 6 quasars reportecP^^ to have warped ac¬ 
cretion disks. Each light curve has been normalized to zero mean and individual curves 
are offset by a constant of 0.5 mag from each other. Error bars shown are standard la 
photometric errors. The CRTS light curve of PG 1302-102 (solid black stars) is also shown 
for comparison. 
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